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Abstract 

We study field strength correlators in presence of a static quark-antiquark pair by 
use of lattice methods. The lattice data have been acquired recently in the context 
of a determination of relativistic corrections to the static interquark potential. We 
extract independent estimates of the form factors of two point correlation func- 
tions and investigate the effect of higher order correlators. Our results confirm the 
dominance of the two point correlators in the regions of intermediate and large 
quark separations, and are compatible with correlation length values obtained from 
different approaches. 
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1 Introduction 



Understanding the QCD dynamics requires insight into non-trivial QCD vac- 
uum properties. In their pioneering work, Shifman, Vainshtein and Zakharov 
[1] focused on imphcations of condensation of quarks and gluons onto ele- 
mentary particle spectroscopy. Soon, it emerged that in general not only the 
condensates had to be considered but also the space-time dependence of cor- 
relation functions as well as a correlation length [2] . Since then, lattice studies 
concentrated on the measurement of the gauge invariant field strength corre- 
lator 

g^{<P{0,x)F^,{x)<P{x,0)Fxp{0)), (1) 



F^, = T»F;„ = TM« , ct>{x, 0) = P eMig J dz,A^(z)}, 

X 

which, at large distances, exhibits an exponential fall-off behaviour that is 
controlled by a correlation length \g. The Stochastic Vacuum Model (SVM) 
offers a systematic way to treat non-perturbative QCD effects which enter 
in terms of gluon field strength correlation functions [3]. In this approach a 
picture emerges in which the correlator (1) plays a basic role in understanding 
colour confinement. It can be decomposed in terms of two form factors D and 
Di, 



g\<P{0,x)F^,{x)<P{x,0)Fxpm 

gHF'm 



{S,x5.p-5^p5,x){D{x^) + D,(x')) 
d 



247VeP(0) + A(0)) 
+{xi^xxS„p - Xp^XpS^x + XuXpSf^x - x^,xxS,^p)—Di{x'^)y (2) 



All parameters within the SVM can in principle be derived in a lattice simu- 
lation directly from the QCD Lagrangian. The simplest linkage between the 
analytic model and lattice gauge theory can be found in the Gaussian approx- 
imation to the SVM. In this approximation cumulants of order higher than 
two in the field strength are neglected in the Wilson loop, such that it is suf- 
ficient to evaluate two-point correlation functions on the lattice. Within this 
approximation, successful predictions ranging from the colour electric field 
distribution within the QCD flux tube between two static quarks and the de- 
termination of the QCD static potential to high energy hadron scattering [4] 
have been produced. To our knowledge all previous lattice evaluations of field 
strength correlators have concentrated on determining two point correlators 
only (see however [5] for an attempt of incorporating higher order effects) 
and relied on Eq. (2) in the simplest choice of the Schwinger fine factors in 
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connection with the coohng method [6]. 

In this paper we have chosen a different approach. We investigate effects of 
field strength correlators in presence of a static qq pair to all orders; these in 
principle contain all physical information about the flux tube configuration 
and the qq heavy quark confined system. With the techniques developed in 
the Wilson loop approach (i.e. manifestly gauge-invariant approach to quark 
dynamics) the complete gauge-invariant quark-antiquark potential at order 
in the quark mass can be obtained from such field strength correlators 
[7,8]. Recently, these analytic expressions were used to obtain the complete 
semirelativistic potential from the lattice [9]. Subsequently, the experimental 
bottomonium spectrum has been fitted to these QCD predictions, and a lattice 
scale as well as the bottom quark mass have been obtained. 

Within the present study, we will use the same (lattice) correlation functions 
1) to investigate the importance of higher order correlation functions within 
the range of interquark distances realized (~ 0.1-1 fm) and 2) to extract in 
this accessible region the behaviour of the D and Di functions. Therefore, 
we do not only aim at a consistency check of the approximation involved 
in the definition of the Gaussian Stochastic Vacuum Model but demonstrate 
indirectly the impact of these non-perturbative correlators onto quarkonia 
spectra. 

The paper is organized as follows. In Sec. 2 we introduce the vacuum expec- 
tation values that have been evaluated on the lattice and their relation to 
the Wilson loop (and therefore to the potential). Furthermore, we outline the 
strategy adopted to identify the contribution of the higher order correlators; 
in Sec. 3 we discuss the numerical results, the validity of the Gaussian approx- 
imation in the SVM and try to extract the large r behaviour of the D and Di 
form factors before we present our conclusions in Sec. 4. 



2 Wilson loop and field strength correlation functions 



Given a set of Schwinger hues $ = {(f){0,u) = cxp{ig j dz^A^{z)}, u G E}, 



connecting points u in Euclidean space E with the origin of the reference 
frame, the Wilson loop (H^(r)) = {exp{ig ^ dz^An{z)}) can be expanded as 

[3,10], 
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iia) r 

n=0 5,(p) 

{Un)(piUn, 0))cum- (3) 

S'(r) denotes a surface with contour T. ^ The cumulants ( )cum are defined in 
terms of path integral averages over gauge fields, 

(0(0, mi)F(mi)0(mi, 0))cum = (0(0, ui)F{ui)<P{u^, 0)), 

(0(0, Mi)F(mi)0(mi, M2)F(m2)0(m2, 0))cum = 

(0(0, «i)F(mi)0(mi, U2)F{U2)<P{U2, 0)) 
- {(f){0,Ui)F{ui)(l){ui,0)) {(f){0,U2)F{u2)(l>{u2,0)), 

where 0(mi,M2) = 0(mi, 0)0(0, ^2) is a shorthand notation. Notice that since 
the left-hand side of Eq. (3) is independent of the set of paths $ as well as 
of the choice of the surface 5", all these dependencies are expected to cancel 
out after summation of the whole scries. However, in general, each cumulant 
appearing within the right-hand side will depend on the choice of 

Let us define 

:F,.x,{x; r) ^ ((/F^.(x)F,,(0)))r - {{gF,^{x)M{gF,,{0)))r. (4) 

This object is required to compute the corrections to the static in- 

terquark potential [7,8]. The double bracket (( ■ ))r stands for the average 
over gauge fields in presence of the Wilson loop Vr(r) 

_ ( ■ WiT)) _ 

//^ {w{r)) ' 

from now on we will assume that the points x = (r, t) and belong to the 
first and second temporal quark line, respectively (see Fig. 1). 

Some particular cases are {Ei = Fi^, = eijkFjk/2): 

£{x; T)^J^uAx; T) = ((E(x) ■ E(0)))r - ((E(,T)))r • ((E(0)))r, 
B{x; r) = -eum€ijk^imjk{x; T) = {{B{x) ■ B(0)))r, 

C{x; r) = ^e,,fc^,4ifc(a;; T) = {{E{x) ■ B(0)))r 

^ Until now it has not been completely proven that the series Eq. (3) has a finite 
convergence radius. We assume that with our choice of $ a region of convergence 
exists. 
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= \tijkT,ulx- r) = ((B(x) ■ E(0)))r, 

Vi{x\ r) = -\^eijkeurnTj/arn{x\ T) = -ej((E(x) A B(0)))r 

= -^tij^eumTimjAix; r) = -ei((B(a;) A E(0)))r, 

where the fact that ((B(x)))r = allows us to disregard some of the discon- 
nected parts. The equalities in the last two relations are due to PT invariance. 
Moreover, due to T invariance, we have Pj(r, t; F) = — Dj(r, —t\ r)[^ 

The function JF can be exactly expressed as a functional derivative of the 
Wilson loop average [7,11]: 



b 5 



J'.uxpix- r) = ,^, iog(iy(r)). (5) 



In Refs. [7,11] the whole semirelativistic qq potential was obtained as a gauge- 
invariant function of the first and second functional derivative on the logarithm 
of the Wilson loop. These derivatives correspond to the vacuum expectation 
values of one or two field strengths insertions in the static Wilson loop [8] that 
have been recently evaluated on the lattice [9]. 

Taking Eq. (3) into account, we obtain, 

J^^uxpix; r) = ^2(0(0, x)F^,{x)(j){x, O)Fap(O)) + n^,xp{x; T), (6) 

where TZ contains all contributions from cumulants of order higher than two: 



n^uxpix; r) = )^ — p ^ j 

n=3 ■ {allperm.(i,j)}5(r) 



X ■ ■ -FapK' = 0) • ■ ■Fp„,.„(M„)0(M„,O))cum. (7) 



In the so-called Stochastic Vacuum Model (SVM) it is assumed that for large 
distances the bilocal cumulant (0(0, x)Fpj^(x)(/)(x, O)Fap(O)) is the dominant 

^ For special geometries of the contour F other general properties can be derived. In 

particular, for a static Wilson loop we obtain {{Eiix))W = = -((^i(0)))r 

r or 

where Vo(r) denotes the static potential between colour sources, separated by a 
spatial distance r. On the lattice this symmetry holds only for large r, such that we 
prefer to determine {{Ei))r by use of the relation ((E(a;)))r • ((E(0)))r = lim ((E(r, t) ■ 

t—>-oo 

E(0)»r. 
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contribution to T in Eq. (6) [4,11,12]. Higher order cumulants are expected 
to be small corrections to the relevant non-perturbative parameters (like the 
string tension) that can be derived from the bilocal one. We have already 
pointed out that this assumption has been successfully tested in potential 
theory and soft high energy scattering (for some reviews see Ref. [4]). However, 
in the present work we aim at a first principles check of this assumption 
on the lattice. The adopted strategy is the following. Lorentz invariance is 
only violated by the contour F within Eq. (6). Therefore, deviations of T 
from Lorentz invariance have to be attributed to the higher order cumulants 
7^ on the right hand side. We interpret the Lorentz invariant part of the 
data as (0(0, x)F^j,(x)0(x, O)-Fap(O)) and parameterize it in terms of two form 
factors. These can subsequently be compared with other lattice determinations 
of these quantities which have been obtained by use of the so-called cooling 
method [6]. We interpret deviations between our data points and continuous 
parameterizations as an estimate for the size of higher order corrections. 



(f, -r/2) 



(f,t + r/2) 




(0,-r/2) (0,0) 



(0,^ + r/2) 



Fig. 1. The rectangular Wilson loop T with two field strength insertions. 



Let us focus on the bilocal cumulant (0(0, x)F^i,{x)(f){x, 0)Fxp{0)) given in Eq. 
(2). For our purposes it is more convenient to define the functions: 



Note that from Eq. (8) the gluon condensate is given by G2 = {g^ (0) / Att'^) — 
QNc/ti"^. The cumulant expansion Eq. (3) takes into account all pertur- 
bative as well as all non-perturbative contributions to the Wilson loop. In 
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particular, for what concerns the bilocal cumulant of Eq. (2), the origin of D 
is purely non-perturbative (with a typical behaviour ~ Q~\^\l^a fQ^- large x) 
while Di contains perturbative (~ l/x"^) as well as non-perturbative contri- 
butions. 

Some special cases of Eq. (2) are, 



/(0(O, x)Fii{x)(t>{x, 0)F,4(0)) = 3(L'x(x') + L>.(x')) 

-2-^D^{x^), (10) 

1.2 

eiimQ,fe/(0(O, x)Fim{x)(t>{x, 0)F,fe(0)) = 12L>x(x') + 8— ^'^(x^), (11) 

ei,fe/(</.(0, x)F,4(x)</.(x, 0)F,fe(0)) = 0, (12) 

e,,fe/(0(O, x)F,fc(x)0(x, 0)F,4(0)) = 0, (13) 

e^,feefe^m5'(0(O,a;)F,•4(a;)0(a;,O)F^^(O)) = -4^D,(a;2), (14) 



ei,feefe,„/((/.(0, x)Fi^{x)(t>{x, 0)F,-4(0)) = -4^L>,(x^). (15) 



With the above definitions, we can rewrite Eq. (6) as 

^nuxp{x] r) = {S^x^^p - Sf,pS^x)D±{x^) 

+ -^{xp,xx5^p - Xf,Xp5^x + x^Xp5^,x - x^xx5p,p)D^{x'^) (16) 
+ T^fiu\p{x;r), 

and, in particular. 



£{x; r) = 3(L'x(x') + D^x^)) - 2^D4x'') + TZuuix; T), (17) 

X 

B{x; r) = 3D^{x^) + 2^D4x^) + \eumeijknimjk{x- T), (18) 

X 4r 

c(x;^) = ^Q,■fe7^^4,•fc(x;^), (19) 
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V{x, r) = -I'j(x; r) = —D^{x ) - -—eijk€kimT^j4im{x; T). (20) 



The lattice data have been produced in order to evaluate the complete 1/m^ 

corrections to the static quark- ant iquark potential [9]. The loop F is the rect- 
angle depicted in Fig. 1 in the limit of large separation in Euclidean time 
between the coordinates of the field strength insertions (ears) and the spa- 
tial closures of the Wilson loop (T » 0, t). In practice, this hmit has been 
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realized by "smearing" the spatial connections within the Wilson loop, i.e. 
by substituting the straight paths by a suitable linear combination of paths 
that maximizes the overlap with the qq ground state (For details see Refs. 
[9,13]). The data have been non-perturbatively renormalized by use of the 
Huntley-Michael method [14]. 

The lattice results provide information on the functions B and V. However, 
not all possible choices of have been realized, such that no data on C is 
available. In the next section we will determine the functions D± and -D*, 
assuming that 71 in Eqs. (17), (18) and (20) can be neglected in the Lorentz- 
symmetric region. Since we have more relations than functions we can also 
provide an independent consistency check of this assumption. In doing so 
we obtain two results. Firstly, as long as the data supply us with space-time 
symmetric functions Z)^ and D^, we can consistently assume that contributions 
from higher order cumulants can be neglected. In fact S, B and T>, depending 
on r, are in general not space-time isotropic functions, whereas D± and D^, are. 
On the other hand, the only source of anisotropy is TZ. Secondly, in the region 
where the bilocal cumulants are established to be the dominant contributions, 
we obtain a parameterization of the form factors that can be compared to 
results already existing in the literature^. 



3 Results 



By inverting these relations we obtain D^, and D± for various values of r and 
t. In order to check space-time symmetry we will consider several combina- 
tions of r and t which yield similar = + t^. The required correlation 
functions S, B and V have been measured on the lattice in the context of the 
extraction of relativistic corrections to the static interquark potential [9]. The 
lattice simulations (with Wilson action, within the valence quark approxima- 
tion to QCD) have been performed at the inverse lattice couplings P = 6.0 and 
f3 = 6.2 on 16^ and 32"^ lattices, respectively. From a fit to the experimental 

^ The definition of the bilocal cumulant is not unique but depends on the set 
of connecting paths <I>. However, we expect the large distance behaviour of the 
cumulants to be universal, in particular the correlation length. 



From Eqs. (17), (18) and (20) we have: 



2rtD^{x^)=x^V{x;T) + ---, 

D^{x') + Id,{x') = i {B{x; F) + S{x; F)) + • • • , 
2 6 

{x^ - At^)D,{x^) = x^ {B{x; F) - S{x; F)) + ■ ■ ■ . 



(21) 
(22) 
(23) 
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bottomonium spectrum, wc obtain inverse lattice spacings of = 2.14 GeV 
and a^^ = 2.94 GeV (a = 0.092 fm and a = 0.067 fm), respectively. In what 
follows, we will neglect systematic errors on these scale estimates from sea 
quark effects and higher order relativistic corrections [9] . 
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0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 0.55 0.6 

x/fm 



Fig. 2. Estimates ofD^ from Eqs. (23) (squares) and (21) (other symbols) at (3 = 6.0. 



Estimates of and from Eqs. (21)-(23) are displayed in Figs. 2 and 
3, respectively. The data points at integer values of t/a in Fig. 2 have been 
obtained from Eq. (23) while the data at half integer values are obtained 
from Eq. (21). We observe that for large x all data collapse onto a universal 
curve. The data from Eq. (23) (squares) provides an independent consistency 
check with respect to Eq. (21). For small x, significant deviations from a 
universal curve are visible in both cases. These deviations signal that higher 
order cumulants become important which is simply an algebraic consequence 
of Eqs. (21) and (23): the function multiplying in Eq. (21) vanishes for 
small values of r = \Jx^ — t"^ and the function multiplying in Eq. (23) 
vanishes at the points x = 2t. For the values of t taken into account, both 
functions are nearly zero for small x. Therefore, in this region, contributions 
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Fig. 3. Estimates of D± from Eqs. (22) and (21) at (3 = 6.0 and P = 6.2 in 
comparison to the asymptotic values (squares). 



of bilocal cumulants to V and 8 — B are suppressed with respect to those 
from higher order cumulants which, as discussed above, are in general not 
space-time symmetric. We find that the observed behaviour is qualitatively 
reproduced (3 — 6.2. In fact, data sets obtained at pairs of similar t-values 
(in physical units) exhibit approximate scaling, such that the data points of 
the figure seem to be relevant to continuum physics, rather than being mere 
lattice artifacts. 

has been obtained by substituting the result on from Eq. (21) (which 
provides us with more precise data than Eq. (23)) into Eq. (22). Contrary to 
the estimates on in this case deviations of data points from an asymptotic 
curve are small, even at small x. In Fig. 3 we display a comparison between 
the t = a data, obtained aX f3 = 6.0 and f3 = 6.2. The results do not only 
weakly depend on t but are almost independent of the lattice spacing a, even 
at these small t values. 
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Within the present statistical errors, universality of the curves is obtained 
from rinin ^ V^a and t^nm ^ 3.5a (tmin ^ 2.5a at P = 6.0) for and 
'"min > \/5a and tmin > 2a for D±{x) onwards. In physical units this distance 
corresponds to Xmin ~ 0.3 fm for and a;niin ~ 0.2 fm for D±, respectively. 




Fig. 4. D± with a one exponential fit to the P = 6.0 data (solid curve) and the curve 
OA/x^ (dashed). 



In Fig. 4 all asymptotic data on D± is displayed, together with a one-expo- 
nential fit to the P — 6.0 data of the form (solid line) 

D^{x')=Aexp{-\x\/Xa), (24) 



The fit parameters remain stable for fit ranges Xmin > 0.3 fm (see Table 1). 
Although even for smaller Xmin, exponential fits yield reasonable values, the 
fit parameters turn out to be unstable under variation of the fit range. Within 
errors the parameters obtained at the two lattice spacings agree with each 
other. We quote as a final result the more precise value obtained a,t P — 6.0, 

A„ = 0.120tn?2fni. (25) 
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Fig. 5. -D* with a one exponential fit to the (3 = 6.0 data. 



Table 1 
Fits to D_L. 

P Xrrnn/a A„/fm A/GeY^ xVNdf 
6.0 3.46 O.UOtil 0.94+11 0.68 
6.2 4.69 0.153t^^ 0.56lf^ 0.80 

In addition, a curve c/x"^, which is the functional form one would expect from 
perturbation theory, with c = 0.4 is plotted in order to show the dominance 
of the perturbative regime in the short range region. 

In Fig. 5 the asymptotic data on is plotted, together with an exponential 
fit of the form 

D,{x'^) = -Bexp{-\x\/Xb), (26) 



to the (3 — 6.0 data. The fit parameters remain stable for fit ranges Xmin > 0.35 
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Table 2 
Fits to D^. 

P Xmin/a Afc/fm ^/GeV^ X^NpF 

6.0 3.77 0.189l^^ 0.47^11 1-25 

6.2 4.92 0.24915^ O.SSIqI 1-47 

fm (see Table 2). Again, within errors the parameters obtained at the two 
lattice agree with each other. We end with the (3 — 6.0 estimate, 

Aft = 0. 1891°:°^^ fm. (27) 

Note that Xa and A;, are lower limits on the corresponding asymptotic correla- 
tion length(s) while A and B are upper limits on the asymptotic amplitudes. 



4 Conclusions 

In this letter we have shown that the lattice data for the long range behaviour 
of the two-field strength insertion on a rectangular Wilson loop arc compatible 
with a two point cumulant approximation. This approximation is crucial in 
the so-called Stochastic Vacuum Model of QCD. Moreover our data are com- 
patible with an exponential fit ~ e"'^'/'*'" with a gluonic correlation length Xg 
of around 0.1-0.2 fm. These findings are in agreement with the quenched re- 
sults obtained in [6] by means of lattice simulations with the cooling method. 
In general one might expect the cooling method, that removes short range 
fluctuations, to have the tendency of underestimating masses and thus to over- 
estimate correlation lengths. On the other hand our results arc more sensitive 
to the perturbative regime which certainly dominates the region in which the 
exponential curves of Fig. 4 and Fig. 5 deviate from the data (typically for 
distances shorter than 0.3 fm). For the moment being, our data suffer under 
much bigger statistical errors than those obtained after cooling. This is natural 
since we do not cut off short range fluctuations. However, we obtain flrst in- 
dependent estimates of the form factors. Unlike cooling, the Hunt ley-Michael 
lattice renormalization does not affect transfer matrix elements. Moreover, it 
has been checked [9] against the exact continuum Gromes [15] and Barchielli- 
Brambilla-Prosperi [8] relations. Another advantage of the present method is 
that the two point cumulant has been extracted from quantities that have a 
direct physical interpretation and play a role in heavy quarkonia spectroscopy. 
We have shown that an evaluation of the bilocal cumulant from existing lattice 
data is possible with the proposed method. We expect that in the near future 
more selective fits can be performed on more precise data. 
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